FAIRE DES CARTES DE FLUX DANS R

L’objectif de cette séance est de montrer comment réaliser des cartes de flux diverses et variées à partir de données matrice pays * pays. Pour cela, le code est souvent assez verbeux. Et beaucoup de choses se font en R base. Il y a beaucoup de manipulation de données. L’idée n’est donc pas ici de commenter tout le code, mais d’expliciter une démarche. C’est à dire montrer comment on peut réaliser des cartes de flux dans R, dans une démarche traçable, partageable et reproductible.

  1. Créez un projet R et un script R.
  2. Créez un repertoire data pour stocker les données.
  3. Créez un répertoire maps dans lequel seront stockées les cartes

Ce document est accessible à l’adresse suivante https://transcarto.github.io/rflows/TRANSCARTO_flows.html

Le code source est disponible ici https://github.com/transcarto/rflows

Les packages

Avant de commencer, voici la liste des packages à installer et à charger. Les 3 packages les plus importants sont sf, mapsf et ttt.

install.packages(sf)
install.packages(remotes)
install.packages(smoothr)
install.packages(readxl)
install.packages(comparator)
install.packages(reshape2)
library(remotes)
install_github("riatelab/mapsf")
install.packages(cartograflow)
install_github("tributetotobler/ttt")
library("sf")
library("mapsf")
library("ttt")
library("readxl")
library("comparator")
library("reshape2")

Les données

Données géométriques

Ici, nous utilisons des données géométriques sur mesure qui permettent de coller exactement avec les données à cartographier. Il s’agit d’un fond de carte des pays du monde dont la nomenclature correspond à celle des données statistiques fournies par les nations unies.

countries <- st_read("https://raw.githubusercontent.com/transcarto/rflows/master/data/world/geom/countries.geojson")
graticule <- st_read("https://raw.githubusercontent.com/transcarto/rflows/master/data/world/geom/graticule.geojson")
bbox <- st_read("https://raw.githubusercontent.com/transcarto/rflows/master/data/world/geom/bbox.geojson")

crs <-
  "+proj=aeqd +lat_0=90 +lon_0=50 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs "
countries <- st_transform(x = countries, crs = crs)
graticule <- st_transform(x = graticule, crs = crs)
bbox <- st_transform(x = bbox, crs = crs)
land <- st_union(countries)

Réalisation d’un template cartographique avec mapsf

col = "#ffc524"
credit = paste0(
  "Françoise Bahoken & Nicolas Lambert, 2021\n",
  "Source: United Nations, Department of Economic\n",
  "and Social Affairs, Population Division (2019)"
)

theme <- mf_theme(
  x = "default",
  bg = "#3b3b3b",
  fg = "#ffc524",
  mar = c(0, 0, 2, 0),
  tab = TRUE,
  pos = "left",
  inner = FALSE,
  line = 2,
  cex = 1.9,
  font = 3
)

template = function(title, file) {
  mf_export(
    countries,
    export = "png",
    width = 1000,
    filename = file,
    res = 96,
    theme = theme,
    expandBB = c(-.02, 0, -.02, 0)
  )
  mf_map(
    bbox,
    col = "#3b3b3b",
    border = NA,
    lwd = 0.5,
    add = TRUE
  )
  mf_map(graticule,
         col = "#FFFFFF50",
         lwd = 0.5,
         add = TRUE)
  mf_map(
    countries,
    col = "#4e4f4f",
    border = "#3b3b3b",
    lwd = 0.5,
    add = TRUE
  )
  # mf_map(links, col = NA,border = "#317691", lwd = 0.5, add = TRUE)
  mf_credits(
    txt = credit,
    pos = "bottomright",
    col = "#1a2640",
    cex = 0.7,
    font = 3,
    bg = "#ffffff30"
  )
  mf_title(title)
}
template("Template cartographique", "maps/template.png")
dev.off()

Données attributaires

Nous utilisons un jeu de données sur les migrations internationales. Migration Stock at subregional level, 2019 Source : United Nations, Department of Economic and Social Affairs, Population Division (2019). Celui-ci est proposé au format xls. Nous l’importons et le mettons en forme via le code ci-dessous.

Voir les données

Téléchargez le fichier UN_MigrantStockByOriginAndDestination_2019.xlsx et placez-le dans votre répertoire data.

Cette opération peut se faire el ligne de code comme suit :

data_url <-
  "https://raw.githubusercontent.com/transcarto/rflows/master/data/world/UN_MigrantStockByOriginAndDestination_2019.xlsx"
file <- "data/UN_MigrantStockByOriginAndDestination_2019.xlsx"
if (!file.exists(file)) {
  download.file(url = data_url, destfile = file)
} 

Choix de la feuille et de l’année de référence

sheet <- "Table 1"
year <- 2019

Import et mise en forme

migr <- data.frame(read_excel(file, skip = 15, sheet = sheet))
migr <- migr[migr[, 1] == year, ]

migr <- migr[!is.na(migr[, 6]), ]
migr <-
  subset(migr,
         select = -c(...1, ...2, ...5, ...4, ...6, Total, Other.North, Other.South))
colnames(migr)[1] <- "i"
migr <- migr[order(migr[, "i"], decreasing = FALSE), ]
for (i in 2:length(colnames(migr))) {
  migr[, i] <- as.numeric(migr[, i])
}

On affecte les codes ISO du fond du carte en ligne et en colonne

ctr <- countries[,2:4] %>% st_drop_geometry()
ctr <- ctr[order(ctr[,"label"], decreasing =FALSE),]
codes <- ctr$adm0_a3_is

# Verification manuelle
ctr$rows <- migr[,"i"]
ctr$cols <- colnames(migr)[-1]
for(i in 1:nrow(ctr)){
  ctr$rows_test[i] = LCS(similarity = TRUE)(ctr$label[i], ctr$rows[i]) / ((nchar(ctr$label[i]) + nchar(ctr$rows[i])) / 2) * 100
  ctr$cols_test[i] = LCS(similarity = TRUE)(ctr$label[i], ctr$cols[i]) / ((nchar(ctr$label[i]) + nchar(ctr$cols[i])) / 2) * 100
}
knitr::kable(ctr[c(0:10),], row.names = F, digits = 1)
un_a3 adm0_a3_is label rows cols rows_test cols_test
4 AFG Afghanistan Afghanistan Afghanistan 100 100.0
8 ALB Albania Albania Albania 100 100.0
12 DZA Algeria Algeria Algeria 100 100.0
16 ASM American Samoa American Samoa American.Samoa 100 92.9
20 AND Andorra Andorra Andorra 100 100.0
24 AGO Angola Angola Angola 100 100.0
660 AIA Anguilla Anguilla Anguilla 100 100.0
28 ATG Antigua and Barbuda Antigua and Barbuda Antigua.and.Barbuda 100 89.5
32 ARG Argentina Argentina Argentina 100 100.0
51 ARM Armenia Armenia Armenia 100 100.0
rownames(migr) <- codes
colnames(migr) <- c("i",codes)
migr <- migr[,-1]
knitr::kable(migr[c(0:15),c(0:15)], row.names = T, digits = 1)
AFG ALB DZA ASM AND AGO AIA ATG ARG ARM ABW AUS AUT AZE BHS
AFG NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
ALB NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
DZA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
ASM NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
AND NA NA NA NA NA NA NA NA 727 NA NA 69 NA NA NA
AGO NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
AIA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
ATG NA NA NA NA NA NA 40 NA NA NA 5 7 NA NA 2
ARG 9 67 105 NA 1 9 NA NA NA 570 1 279 1039 NA NA
ARM NA NA NA NA NA NA NA NA NA NA NA NA NA 78478 NA
ABW NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
AUS 59798 3872 1745 0 31 592 NA 63 17583 1629 74 NA 18091 740 327
AUT 20561 3715 1522 NA 2 402 NA 8 1797 3601 NA 2939 NA 1270 25
AZE 178 NA NA NA NA NA NA NA NA 142650 NA NA NA NA NA
BHS NA NA NA NA NA NA NA 14 117 NA NA 88 42 NA NA

On transpose la matrice

migr <- t(migr)

Avec reshape2, on la convertit au format i,j,fij

migr <- melt(migr)
colnames(migr) = c("i","j","fij")
migr <- migr[!is.na(migr$fij),]
migr = migr[migr$fij>0,]
migr = migr[order(migr$fij, decreasing = TRUE),]
knitr::kable(migr[c(0:10),], row.names = F, digits = 1)
i j fij
MEX USA 11489684
SYR TUR 3743494
IND ARE 3419875
RUS UKR 3308515
UKR RUS 3269248
BGD IND 3103664
CHN USA 2899267
IND USA 2661470
KAZ RUS 2559711
RUS KAZ 2458414

Sauvegarder le fichier mis en forme au format csv

write.csv(migr, "data/migr.csv", row.names = FALSE)

Ce fichier de données correctement formaté est dorénavant accesible comme ceci.

migr <- read.csv("data/migr.csv")

Tout est prêt. Avançons…

Premières explorations

L’effet Spaghetti

links <-
  mf_get_links(
    x = countries,
    df = migr,
    x_id = "adm0_a3_is",
    df_id = c("i", "j")
  )
template("L'effet Spaghetti ", "maps/spaghetti.png")
mf_map(links, col = col, add = TRUE)
mf_map(land,
       col = NA,
       border = "#3b3b3b",
       add = TRUE)
dev.off()

Un pays de référence

Pour simplifier la carte, on peut choisir un seuk pays de référence

ISO3 <- "FRA"
label = "France"

Jointure et mise en forme des données

countr <- countries[, c("adm0_a3_is", "label")]
migrFRA <- migr[migr$j == ISO3, ]
migrFRA$fij <- as.numeric(migrFRA$fij)
maxval = max(migrFRA$fij)
total = round(sum(migrFRA$fij) / 1000000,1)
countr <-
  merge(
    x = countr,
    y = migrFRA,
    by.x = "adm0_a3_is",
    by.y = "i",
    all.x = TRUE
  )
countr <- countr[-3]
colnames(countr) <- c("id", "label", "fij", "geometry")
knitr::kable(countr[c(0:10),], row.names = F, digits = 1)
id label fij geometry
ABW Aruba 11 MULTIPOLYGON (((-7476945 42…
AFG Afghanistan 6887 MULTIPOLYGON (((2474775 -53…
AGO Angola 23438 MULTIPOLYGON (((-4917506 -1…
AIA Anguilla 10 MULTIPOLYGON (((-7351488 31…
ALB Albania 7371 MULTIPOLYGON (((-2639654 -4…
AND Andorra 1079 MULTIPOLYGON (((-3952645 -3…
ARE United Arab Emirates 862 MULTIPOLYGON (((785851 -712…
ARG Argentina 14253 MULTIPOLYGON (((-14113355 7…
ARM Armenia 21012 MULTIPOLYGON (((-348529.3 -…
ASM American Samoa 1 MULTIPOLYGON (((7561304 878…

Une première carte simple

template(paste0("En 2019, il y avait ",total, " millions d'étrangers en France"),
         "maps/prop1.png")
#mf_map(countr[countr$id == ISO3,], col = col, border = "red", lwd = 2, add = TRUE)
mf_map(
  countr[countr$id != ISO3, ],
  var = "fij",
  col = col,
  border = "white",
  type = "prop",
  val_max = maxval,
  inches = 0.4,
  leg_title_cex = 1.2,
  leg_val_cex   = 0.8,
  leg_pos = "bottomleft",
  leg_title = "Nombre de personnes"
)
mf_map(
  countr[countr$id == ISO3, ],
  col = NA,
  border = "#e36019",
  lwd = 2,
  add = TRUE
)
dev.off()

La carte symétrique

countr <- countries[, c("adm0_a3_is", "label")]
migrFRA <- migr[migr$i == ISO3, ]
migrFRA$fij <- as.numeric(migrFRA$fij)
total = round(sum(migrFRA$fij) / 1000000,1)
countr <-
  merge(
    x = countr,
    y = migrFRA,
    by.x = "adm0_a3_is",
    by.y = "j",
    all.x = TRUE
  )
countr <- countr[-3]
colnames(countr) <- c("id", "label", "fij", "geometry")
template(paste0("En 2019, il y avait ",total, " millions de Français à l'étranger"),
         "maps/prop2.png")
mf_map(
  countr[countr$id != ISO3, ],
  var = "fij",
  col = col,
  border = "white",
  type = "prop",
  val_max = maxval,
  inches = 0.4,
  leg_title_cex = 1.2,
  leg_val_cex   = 0.8,
  leg_pos = "bottomleft",
  leg_title = "Nombre de personnes"
)
mf_map(
  countr[countr$id == ISO3, ],
  col = NA,
  border = "#e36019",
  lwd = 2,
  add = TRUE
)
dev.off()

On peut faire la même carte en faisant varier l’épaisseur des liens

ISO3 <- "FRA"
label = "France"
migrtoFRA <- migr[migr$j == ISO3,]
migrtoFRA$fij <- as.numeric(migrtoFRA$fij)
links <-
  mf_get_links(
    x = countries,
    df = migrtoFRA,
    x_id = "adm0_a3_is",
    df_id = c("i", "j")
  )
template(
  paste0("Origine des personnes migrantes vivant en ", label, " en 2019"),
  "maps/links1.png"
)
mf_map(
  links,
  var = "fij",
  col = col,
  border = "white",
  type = "prop",
  inches = 10,
  leg_title_cex = 1.2,
  leg_val_cex   = 0.8,
  leg_pos = "bottomleft",
  leg_title = "Nombre de personnes"
)
mf_map(
  countries[countries$adm0_a3_is == ISO3,],
  col = "#4e4f4f",
  border = col,
  lwd = 1.5,
  add = TRUE
)
dev.off()

Filtrages et indicateurs

FRANCOISE, TU METS TES TRUCS ICI !

Lorem ipsum dolor sit amet, consectetur adipiscing elit. Vivamus ac commodo ante. Sed tincidunt tincidunt sollicitudin. Mauris odio orci, viverra et porttitor vitae, porta et arcu. Nam quis neque at elit accumsan fringilla. Cras rhoncus efficitur malesuada. Donec auctor, mauris sit amet mollis dignissim, purus urna ultricies ligula, non dictum sem quam in justo. Maecenas sit amet est accumsan, ullamcorper nulla luctus, efficitur diam. Aenean ac magna ut enim lacinia ultrices nec et sapien. Praesent scelerisque massa eros, vel tempus orci maximus at. Duis eget ipsum auctor, luctus ante in, egestas libero. Vestibulum vehicula ex a aliquam aliquam. Donec non efficitur risus. Aenean ut venenatis nisi, vitae iaculis nibh. Cras fermentum orci vel tempor sodales.

Sed viverra ut ipsum in commodo. Quisque tempus tempus tortor ut feugiat. Morbi viverra, metus id feugiat vulputate, augue libero condimentum ex, vestibulum volutpat odio metus id arcu. Integer ullamcorper sed sapien ut sollicitudin. Ut aliquet, leo at elementum ornare, sem ante lacinia ante, ac porttitor justo arcu nec dolor. Cras porta nisl lobortis leo dignissim ullamcorper. Cras eu lorem imperdiet, malesuada risus ac, tempus neque. Integer in erat consequat, posuere sapien quis, pulvinar turpis. In sagittis cursus commodo. Praesent pellentesque commodo velit, quis suscipit dui sodales tristique. Mauris eleifend quam et odio viverra, quis suscipit ex semper. In semper id sapien id egestas. Mauris eros metus, rhoncus eu arcu convallis, sagittis tempor nunc. Pellentesque hendrerit, tortor at lacinia lacinia, neque neque interdum lacus, in finibus tortor metus et arcu. Praesent sed viverra lectus, nec elementum velit. In faucibus neque in risus ultricies cursus vel sed magna.

Mauris a ante nec mi ornare egestas sit amet vel mauris. Maecenas ac dolor id dolor facilisis fermentum id a orci. Praesent sed dolor non nisl vulputate pulvinar. Donec vehicula vitae massa vel semper. Sed sit amet cursus odio. Fusce blandit ligula mollis justo consectetur, eget finibus nulla molestie. Morbi convallis nulla non mi finibus tempor. Vestibulum sagittis vitae mauris ut pulvinar. Morbi aliquam iaculis leo. Cras massa odio, commodo eu libero sit amet, dictum condimentum dolor. Donec posuere rutrum purus vitae euismod. Phasellus vel leo nec nisl varius luctus. Cras sed suscipit quam.

Vers des cartes un peu plus graphiques

Une carte un peu plus sophistiquée avec packcircles

Avec le code ci-dessous, on cherche à réaliser une carte à la façon de cette application interactive.

https://analytics.huma-num.fr/Nicolas.Lambert/migrexplorer/

ISO3 <- "FRA"
label = "France"
migrFRA <- migr[migr$j == ISO3,]
migrFRA$fij <- as.numeric(migrFRA$fij)
migrFRA <-
  rbind.data.frame(migrFRA, c(
    i = ISO3,
    j = ISO3,
    fij = sum(as.numeric(migrFRA$fij))
  ))
countr <- countries[, "adm0_a3_is"]
countr <-
  merge(
    x = countr,
    y = migrFRA,
    by.x = "adm0_a3_is",
    by.y = "i",
    all.x = TRUE
  )
colnames(countr) <- c("i", "j", "fij", "geometry")
knitr::kable(countr[c(0:10),], row.names = F, digits = 1)
i j fij geometry
ABW FRA 11 MULTIPOLYGON (((-7476945 42…
AFG FRA 6887 MULTIPOLYGON (((2474775 -53…
AGO FRA 23438 MULTIPOLYGON (((-4917506 -1…
AIA FRA 10 MULTIPOLYGON (((-7351488 31…
ALB FRA 7371 MULTIPOLYGON (((-2639654 -4…
AND FRA 1079 MULTIPOLYGON (((-3952645 -3…
ARE FRA 862 MULTIPOLYGON (((785851 -712…
ARG FRA 14253 MULTIPOLYGON (((-14113355 7…
ARM FRA 21012 MULTIPOLYGON (((-348529.3 -…
ASM FRA 1 MULTIPOLYGON (((7561304 878…

Cercles avec packcircles (Dorling style)

library(packcircles)
dots = countr
st_geometry(dots) <-
  st_centroid(sf::st_geometry(dots), of_largest_polygon = TRUE)
dots <- data.frame(dots$i, dots["fij"], st_coordinates(dots))
dots = dots[, c("dots.i", "X", "Y", "fij")]
colnames(dots) <- c("id", "x", "y", "v")
dots <- dots[!is.na(dots$v), ]

k = 700000 # pour ajuster la taille des cercles
itermax = 10 # nombre d'iterations
delta = 35000
dat.init <- dots[, c("x", "y", "v", "id")]
dat.init$v <- sqrt(as.numeric(dat.init$v) * k)
simulation <- circleRepelLayout(
  x = dat.init,
  xysizecols = 1:3,
  wrap = FALSE,
  sizetype = "radius",
  maxiter = itermax,
  weights = 1
)$layout
circles <- st_buffer(sf::st_as_sf(
  simulation,
  coords = c('x', 'y'),
  crs = sf::st_crs(countries)
),
dist = simulation$radius - delta)

circles$v = dots$v
circles$id = dots$id

Links

# Links

dots$j = "FRA"

links <-
  mf_get_links(
    x = circles,
    df = migrFRA,
    x_id = "id",
    df_id = c("i", "j")
  )
links$fij = as.numeric(links$fij)

Réalisation de la carte

template("Les étrangers en France, 2019", "maps/migrexplorer1.png")

col2 = "#4e4f4f"

mf_map(
  land,
  col = "#4e4f4f",
  border = "#3b3b3b",
  lwd = 0.5,
  add = TRUE
)


mf_map(
  links,
  var = "fij",
  col = col,
  border = "#3b3b3b",
  type = "prop",
  lwd_max = 160,
  leg_pos = "n",
  add = TRUE
)

mf_map(
  circles[circles$id != ISO3, ],
  var = "fij",
  col = col,
  border = "#3b3b3b",
  lwd = 1.5,
  add = TRUE
)

mf_map(
  circles[circles$id == ISO3, ],
  var = "fij",
  col = col2,
  border = col,
  lwd = 2.5,
  add = TRUE
)

t =  circles[circles$id != ISO3, ]
mf_label(
  t,
  var = "id",
  halo = FALSE,
  cex = sqrt(as.numeric(t$v) / 1200000),
  col = col2,
  overlap = TRUE,
  lines = FALSE
)

t =  circles[circles$id == ISO3, ]
mf_label(
  t,
  var = "id",
  halo = FALSE,
  cex = sqrt(as.numeric(t$v) / 1200000),
  col = col,
  overlap = TRUE,
  lines = FALSE
)

dev.off()

Comme précédemment, on peut faire la carte en symétrie en inversant i et j.

ISO3 <- "FRA"
label = "France"
migrFRA <- migr[migr$i == ISO3,] # ici
migrFRA$fij <- as.numeric(migrFRA$fij)
migrFRA <-
  rbind.data.frame(migrFRA, c(
    i = ISO3,
    j = ISO3,
    fij = sum(as.numeric(migrFRA$fij))
  ))
countr <- countries[, "adm0_a3_is"]
countr <-
  merge(
    x = countr,
    y = migrFRA,
    by.x = "adm0_a3_is",
    by.y = "j", # là
    all.x = TRUE
  )
colnames(countr) <- c("i", "j", "fij", "geometry")
dots = countr
st_geometry(dots) <-
  st_centroid(sf::st_geometry(dots), of_largest_polygon = TRUE)
dots <- data.frame(dots$i, dots["fij"], st_coordinates(dots))
dots = dots[, c("dots.i", "X", "Y", "fij")]
colnames(dots) <- c("id", "x", "y", "v")
dots <- dots[!is.na(dots$v), ]

k = 700000 # pour ajuster la taille des cercles
itermax = 10 # nombre d'iterations
delta = 35000
dat.init <- dots[, c("x", "y", "v", "id")]
dat.init$v <- sqrt(as.numeric(dat.init$v) * k)
simulation <- circleRepelLayout(
  x = dat.init,
  xysizecols = 1:3,
  wrap = FALSE,
  sizetype = "radius",
  maxiter = itermax,
  weights = 1
)$layout
circles <- st_buffer(sf::st_as_sf(
  simulation,
  coords = c('x', 'y'),
  crs = sf::st_crs(countries)
),
dist = simulation$radius - delta)

circles$v = dots$v
circles$id = dots$id

Links

# Links

dots$j = "FRA"

links <-
  mf_get_links(
    x = circles,
    df = migrFRA,
    x_id = "id",
    df_id = c("i", "j")
  )
links$fij = as.numeric(links$fij)

Réalisation de la carte

template("Les français à l'étranger, 2019", "maps/migrexplorer2.png")

col2 = "#4e4f4f"

mf_map(
  land,
  col = "#4e4f4f",
  border = "#3b3b3b",
  lwd = 0.5,
  add = TRUE
)


mf_map(
  links,
  var = "fij",
  col = col,
  border = "#3b3b3b",
  type = "prop",
  lwd_max = 160,
  leg_pos = "n",
  add = TRUE
)

mf_map(
  circles[circles$id != ISO3, ],
  var = "fij",
  col = col,
  border = "#3b3b3b",
  lwd = 1.5,
  add = TRUE
)

mf_map(
  circles[circles$id == ISO3, ],
  var = "fij",
  col = col2,
  border = col,
  lwd = 2.5,
  add = TRUE
)

t =  circles[circles$id != ISO3, ]
mf_label(
  t,
  var = "id",
  halo = FALSE,
  cex = sqrt(as.numeric(t$v) / 1200000),
  col = col2,
  overlap = TRUE,
  lines = FALSE
)

t =  circles[circles$id == ISO3, ]
mf_label(
  t,
  var = "id",
  halo = FALSE,
  cex = sqrt(as.numeric(t$v) / 1200000),
  col = col,
  overlap = TRUE,
  lines = FALSE
)

dev.off()

Ces cartes, on peut les retrouver dans l’application MigrExplorer mise en ligne via R shiny.

https://gitlab.huma-num.fr/nlambert/migrexplorer/-/tree/master

Changer de maillage

Contrairement aux cartes pays * pays, cartographier les flux au niveau régional permet de mieux percevoir la logique des mobilités internationales. Cette carte, pas très élégantes, a été réalisée et présentée par François Héron pour ses cours au Collège de France.

Et si on esseyait de la reproduire en R ?

Pour celà, nous fabriquons des données au niveau subrégional à partir d’une clé d’aggrégations contenu dans le ficher countries.

knitr::kable(countries[c(0:10),c("adm0_a3_is", "label","Code2","Label2")], row.names = F, digits = 1)
adm0_a3_is label Code2 Label2 geometry
BGR Bulgaria 923 Eastern Europe MULTIPOLYGON (((-1882818 -4…
MMR Myanmar 920 South-Eastern Asia MULTIPOLYGON (((5416951 -56…
BDI Burundi 910 Eastern Africa MULTIPOLYGON (((-3418256 -9…
BLR Belarus 923 Eastern Europe MULTIPOLYGON (((-1406024 -3…
KHM Cambodia 920 South-Eastern Asia MULTIPOLYGON (((7198820 -51…
DZA Algeria 912 Northern Africa MULTIPOLYGON (((-3911770 -4…
CMR Cameroon 911 Middle Africa MULTIPOLYGON (((-5196562 -7…
CAN Canada 918 Northern America MULTIPOLYGON (((-2925928 15…
CPV Cabo Verde 914 Western Africa MULTIPOLYGON (((-7996256 -2…
CYM Cayman Islands 915 Caribbean MULTIPOLYGON (((-5899896 51…

Géométries

subregions <-
  aggregate(countries, by = list(countries$Code2), FUN = head, 1)
subregions <- subregions[, c("Code2", "Label2")]
st_geometry(subregions) <-
  st_cast(subregions$geometry, "MULTIPOLYGON")
colnames(subregions) <- c("id", "label", "geometry")
template("Subregions", "maps/subregions.png")
mf_map(
  subregions,
  col = "#4e4f4f",
  border = col,
  lwd = 0.5,
  add = TRUE
)
mf_label(
  x = subregions,
  var = "label",
  halo = TRUE,
  bg = "#4e4f4f",
  cex = 0.8,
  col = col,
  overlap = TRUE,
  lines = FALSE
)
dev.off()

Données attributaires

keys <- data.frame(countries[, c("adm0_a3_is", "Code2")])
keys$geometry <- NULL
migr <- merge(x = migr,
              y = keys,
              by.x = "i",
              by.y = "adm0_a3_is")
colnames(migr)[4] <- "subreg_i"
migr <- merge(x = migr,
              y = keys,
              by.x = "j",
              by.y = "adm0_a3_is")
colnames(migr)[5] <- "subreg_j"
migr$id <- paste0(migr$subreg_i, "_", migr$subreg_j)
migr2 <- aggregate(migr$fij, by = list(migr$id), FUN = sum)
migr2$i <- sapply(strsplit(migr2$Group.1, "_"), "[", 1)
migr2$j <- sapply(strsplit(migr2$Group.1, "_"), "[", 2)
migr2 <- migr2[, c("i", "j", "x")]
colnames(migr2)[3] <- "fij"
migr2$fij <- round(migr2$fij / 1000, 0)
knitr::kable(migr2[c(0:10),], row.names = F, digits = 1)
i j fij
5500 5500 483
5500 5501 12
5500 906 28
5500 912 4
5500 913 0
5500 914 2
5500 915 0
5500 916 0
5500 918 137
5500 922 95

On ajoute au fond de carte les flux intrarégionaux

flowsintra <- migr2[migr2$i == migr2$j,c("i","fij")]
colnames(flowsintra) <- c("id","intra")
subregions <- merge(x = subregions, y = flowsintra, by = "id")
knitr::kable(subregions[c(0:10),], row.names = F, digits = 1)
id label intra geometry
906 Eastern Asia 5202 MULTIPOLYGON (((6897029 -39…
910 Eastern Africa 5330 MULTIPOLYGON (((1139915 -12…
911 Middle Africa 1537 MULTIPOLYGON (((-4988352 -1…
912 Northern Africa 351 MULTIPOLYGON (((-1695160 -7…
913 Southern Africa 715 MULTIPOLYGON (((-3980523 -1…
914 Western Africa 6625 MULTIPOLYGON (((-5014884 -7…
915 Caribbean 864 MULTIPOLYGON (((-8056812 28…
916 Central America 641 MULTIPOLYGON (((-7273542 55…
918 Northern America 1114 MULTIPOLYGON (((-1560126 -6…
920 South-Eastern Asia 6856 MULTIPOLYGON (((6961013 -67…

Calcul des interactions inter régionales (A -> B) + (B -> A)

migr2 <- migr2[migr2$i != migr2$j,]
for (k in 1:length(migr2$i)) {
  val1 <- migr2$fij[k]
  val2 <-
    migr2[migr2$i == migr2$j[k] & migr2$j == migr2$i[k], "fij"]
  migr2$interaction[k] <- sum(val1, val2)
}

# Suppression des doublons
interactions = data.frame(matrix(
  ncol = 3,
  nrow = 0,
  dimnames = list(NULL, c("i", "j", "interaction"))
))
for (k in 1:length(migr2$i)) {
  idi = migr2$i[k]
  idj = migr2$j[k]
  test = length(interactions[(interactions$i == idi &
                                interactions$j == idj) |
                               (interactions$i == idj & interactions$j == idi), "interaction"])
  if (test == 0) {
    interactions <-
      rbind(interactions, data.frame(
        i = idi,
        j = idj,
        interaction = migr2$interaction[k]
      ))
  }
}
knitr::kable(interactions[c(0:10),], row.names = F, digits = 1)
i j interaction
5500 5501 28
5500 906 130
5500 912 4
5500 913 0
5500 914 2
5500 915 0
5500 916 0
5500 918 137
5500 922 261
5500 923 9999

On élimine les petits flux

threshold <- 2000
interactions <- interactions[interactions$interaction >= threshold,]

Calcul des liens

links <-
  mf_get_links(
    x = subregions,
    df = interactions,
    x_id = "id",
    df_id = c("i", "j")
  )

Cartographie

template("L'Arique, un continent encore isolé dans la mondialisation", "maps/heran.png")

col2 = "#4e4f4f"

mf_map(
  subregions,
  col = "#4e4f4f",
  border = "#3b3b3b",
  lwd = 0.5,
  add = TRUE
)


mf_map(
  links,
  var = "interaction",
  col = col,
  border = "#3b3b3b",
  type = "prop",
  lwd_max = 25,
  leg_pos = "bottomleft",
  leg_title = paste0("Migratons INTER régionales (interactions)\n(A -> B) + (B -> A)\nSeuil : ",threshold, "\nen milliers de personnes"),
  add = TRUE
)

mf_map(
  subregions,
  var = "intra",
  col = "#3b3b3b",
  border = col,
  lwd = 1.5,
  type = "prop",
  symbol = "square",
  leg_pos = "topright",
  leg_title = "Migrations INTRA\nrégionale nen 2019\n(en milliers)",
  add = TRUE
)

mf_label(
  subregions,
  var = "intra",
  halo = FALSE,
  cex = sqrt(as.numeric(subregions$intra) / 12000),
  col = col,
  overlap = TRUE,
  lines = FALSE
)

mf_label(
  links,
  var = "interaction",
  halo = TRUE,
  cex = 0.5,
  col = col2,
  bg = col,
  r = 0.1,
  overlap = FALSE,
  lines = FALSE
)

dev.off()

Problème : avec seulement mapsf, on a du mal à représenter des flêches et surtout, à la fois des flêches A -> B et B -> A. La solution : Flowmapper 👍

Flowmapper

flowmapper() est une fonction du package ttt (en cours de développement).

library(ttt)

Les données

Dans le package ttt, il y a des données d’exemple au niveau subrégional. Chargeons-les.

subregions <- st_read(system.file("subregions.gpkg", package="flowmapper")) %>% st_transform(crs)
migr <- read.csv(system.file("migrantstocks2019.csv", package="flowmapper"))

On ne consrve que les flux importants

threshold <- 1500
migr <- migr[migr$fij >= threshold, ]
knitr::kable(migr[c(0:10),], row.names = F, digits = 1)
i j fij
5500 923 5603
5501 5501 11177
5501 918 5334
5501 920 1666
5501 922 18402
5501 924 2551
906 906 5202
906 918 5700
910 910 5330
910 913 1538
flows <- ttt_flowmapper(
  x = subregions,
  xid = "id",
  df = migr,
  dfid = c("i", "j"),
  dfvar = "fij",
  plot = FALSE
)

Liens

template("ttt_flowmapper$links", "maps/ttt_links.png")
mf_map(
  subregions,
  col = "#4e4f4f",
  border = "#3b3b3b",
  lwd = 0.5,
  add = TRUE
)
mf_map(flows$links,
       col = col,
       lwd = 3,
       add = TRUE)
dev.off()

Cercles

template("ttt_flowmapper$circles", "maps/ttt_circles.png")
mf_map(
  subregions,
  col = "#4e4f4f",
  border = "#3b3b3b",
  lwd = 0.5,
  add = TRUE
)
mf_map(flows$circles, col = col, add = TRUE)
dev.off()

Flêches

template("ttt_flowmapper$flows", "maps/ttt_flows.png")
mf_map(
  subregions,
  col = "#4e4f4f",
  border = "#3b3b3b",
  lwd = 0.5,
  add = TRUE
)
mf_map(flows$flows, col = col, add = TRUE)
dev.off()

Visualisation par défaut

template("flowmappze", "maps/ttt_flows.png")
mf_map(
  subregions,
  col = "#4e4f4f",
  border = "#3b3b3b",
  lwd = 0.5,
  add = TRUE
)
flows <- ttt_flowmapper(
  x = subregions,
  xid = "id",
  type = "arrows",
  df = migr,
  dfid = c("i", "j"),
  dfvar = "fij",
  col = col,
  border = "#424242",
  border2 = col,
  add = TRUE
)
ttt_flowmapperlegend(x = flows, title = "Flux", col = col)
dev.off()

La VV taille, c’est aussi la surface

template("La surface des fleches", "maps/ttt_surface.png")
mf_map(
  subregions,
  col = "#4e4f4f",
  border = "#3b3b3b",
  lwd = 0.5,
  add = TRUE
)
ttt_flowmapper(
  x = subregions,
  xid = "id",
  type = "arrows",
  size = "area",
  df = migr,
  dfid = c("i", "j"),
  dfvar = "fij",
  col = col,
  border = "#424242",
  border2 = col,
  add = TRUE
)
dev.off()

Epaisseur vs Surface

Interactions (type = “rect”)

migr2 <- data.frame(i = integer(), j = integer(), fij = integer())

for (k in 1:length(migr$i)) {
  val1 <- migr$fij[k]
  val2 <- migr[migr$i == migr$j[k] & migr$j == migr$i[k], "fij"]
  val <- sum(val1, val2)
  idi =  migr$i[k]
  idj =  migr$j[k]
  test <-
    length(migr2[(migr2$i == idi &
                    migr2$j == idj) | (migr2$i == idj & migr2$j == idi), "fij"])
  if (test == 0) {
    migr2 <- rbind(migr2, data.frame(i = idi, j = idj, fij = val))
  }
}
migr2 <- migr2[migr2$i != migr2$j, ] 
head(migr2)
##      i   j   fij
## 1 5500 923  9999
## 3 5501 918  5334
## 4 5501 920  3221
## 5 5501 922 18402
## 6 5501 924  2551
## 8  906 918  5700
template("tInteractions", "maps/ttt_interactions.png")
c <- ttt_flowmapper(
  x = subregions,
  xid = "id",
  size = "thickness",
  type = "rect",
  df = migr2,
  dfid = c("i", "j"),
  dfvar = "fij",
  col = col,
  border = "#424242",
  border2 = col,
  add = TRUE
)
dev.off()

Combiner flux intra et flux inter

intra <- migr[migr$i == migr$j, ]
intra <- intra[, c("i", "fij")]
colnames(intra) <- c("id", "nb")
knitr::kable(intra, row.names = F, digits = 1)

Calcul des flux (plot = FALSE)

flows <- ttt_flowmapper(
  x = subregions,
  xid = "id",
  df = migr,
  dfid = c("i", "j"),
  dfvar = "fij",
  size = "thickness",
  type = "arrows",
  decreasing = FALSE,
  add = TRUE,
  lwd = 1,
  col = col,
  border = "#424242",
  k = NULL,
  k2 = 60,
  df2 = intra,
  df2id = "id",
  df2var = "nb",
  col2 = "#3b3b3b",
  border2 = col,
  plot = FALSE
)

Affichage de la carte avec mapsf

template("Flux inter et flux intra", "maps/interintra.png")

mf_shadow(x = flows$flows, col = "grey70", cex = 0.2, add = TRUE)

mf_map(
  flows$flows,
  var = "fij",
  col = col,
  border = "#3b3b3b",
  leg_pos = "n",
  add = TRUE
)

mf_map(
  flows$circles,
  var = "fij",
  col = "#3b3b3b",
  border = col,
  lwd = 1.5,
  leg_pos = "n",
  add = TRUE
)

mf_label(
  flows$circles,
  var = "nb",
  halo = FALSE,
  cex = sqrt(as.numeric(flows$circles$nb) / 18000),
  #cex = 1,
  col = col,
  overlap = TRUE,
  lines = FALSE
)


mf_label(
  flows$flows,
  var = "fij",
  halo = TRUE,
  cex = 0.7,
  col = col2,
  bg = col,
  r = 0.1,
  overlap = FALSE,
  lines = FALSE
)



dev.off()

Reprojection

1 - calcul en projection polaire

tmp <- ttt_flowmapper(
  x = subregions,
  xid = "id",
  type = "arrows",
  df = migr,
  dfid = c("i", "j"),
  dfvar = "fij",
  col = "#ffc524",
  border = "#424242",
  border2 = "#ffc524",
  plot = FALSE
)

2 - reprojection & nouveau template

crs <-
  "+proj=ortho +lat_0=42.5333333333 +lon_0=-72.53333333339999 +x_0=0 +y_0=0 +a=6370997 +b=6370997 +units=m +no_defs"
flows <- smoothr::densify(tmp$flows, n = 30) %>% st_transform(crs)
dots <- st_transform(tmp$circles, crs)
subregions <- st_transform(subregions, crs)
graticule <- st_transform(graticule, crs)
bbox <- st_transform(bbox, crs)

3 - affichage

title = "Flux sur Globe"
file =   "maps/ttt_globe.png"

mf_export(
  subregions,
  export = "png",
  width = 1000,
  filename = file,
  res = 96,
  theme = theme,
  expandBB = c(-.02, 0,-.02, 0)
)

mf_map(
  bbox,
  col = "#3b3b3b",
  border = NA,
  lwd = 0.5,
  add = TRUE
)

mf_map(graticule,
       col = "#FFFFFF50",
       lwd = 0.5,
       add = TRUE)

mf_map(
  subregions,
  col = "#4e4f4f",
  border = "#3b3b3b",
  lwd = 0.5,
  add = TRUE
)

mf_credits(
  txt = credit,
  pos = "bottomright",
  col = "#1a2640",
  cex = 0.7,
  font = 3,
  bg = "#ffffff30"
)

mf_map(flows, col = col, add = TRUE)

mf_map(dots, col = col, add = TRUE)

mf_title(title)

dev.off()

Vers une implémantation dans Observable (svg/d3js)

Ouvrir l’application

A vous de jouer

Et si on essayait de faire des cartes de flux sur un fond de carte déformé. Ici, par la population en 2019.

Les données

migrCountries <- read.csv("https://raw.githubusercontent.com/transcarto/rflows/master/data/world/cartogram/migr2019_T.csv")
migrSubregions <-  read.csv("https://raw.githubusercontent.com/transcarto/rflows/master/data/world/cartogram/migrSubregions2019_T.csv")

Les géométries

countriesPop <- st_read("https://raw.githubusercontent.com/transcarto/rflows/master/data/world/cartogram/countriesPop.geojson")
subregionsPop <- st_read("https://raw.githubusercontent.com/transcarto/rflows/master/data/world/cartogram/regionsPop.geojson")
gridPop <- st_read("https://raw.githubusercontent.com/transcarto/rflows/master/data/world/cartogram/grid.geojson")

Quelques variables d’affichage

col = "#ffc524"
votrenom = "Miles Davis & Frances Taylor, kings of cool, 1965"
credit = paste0(
  votrenom,"\n",
  "Source: United Nations, Department of Economic\n",
  "and Social Affairs, Population Division (2019)"
)

theme <- mf_theme(
  x = "default",
  bg = "#3b3b3b",
  fg = "#ffc524",
  mar = c(0, 0, 2, 0),
  tab = TRUE,
  pos = "left",
  inner = FALSE,
  line = 2,
  cex = 1.9,
  font = 3
)

template = function(title, file) {
  mf_export(
    countriesPop,
    export = "png",
    width = 1000,
    filename = file,
    res = 96,
    theme = theme,
    expandBB = c(-.02, 0, -.02, 0)
  )

  mf_map(gridPop,
         col = "#FFFFFF70",
         lwd = 0.4,
         add = TRUE)
  mf_map(
    countriesPop,
    col = "#4e4f4f",
    border = "#3b3b3b",
    lwd = 0.5,
    add = TRUE
  )
  
    mf_map(
    subregionsPop,
    col = "NA",
    border = col,
    lwd = 0.5,
    add = TRUE
  )
  # mf_map(links, col = NA,border = "#317691", lwd = 0.5, add = TRUE)
  mf_credits(
    txt = credit,
    pos = "bottomright",
    col = "#1a2640",
    cex = 0.5,
    font = 3,
    bg = "#ffffff30"
  )
  mf_title(title)
}

A vous de jouer…

template("World Population, 2019", "maps/cartogram.png")

# METTEZ DES TRUCS ICI ! 

dev.off()